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Abstract 

We present thermodynamic relationships between the free energy of the phase-field crystal (PFC) 
model and thermodynamic state variables for bulk phases under hydrostatic pressure. This rela¬ 
tionship is derived based on the thermodynamic formalism for crystalline solids of Larche and 
Calm [Larche and Cahn, Acta Metallurgica, Vol. 21, 1051 (1973)]. We apply the relationship to 
examine the thermodynamic processes associated with varying the input parameters of the PFC 
model: temperature, lattice spacing, and the average value of the PFC order parameter, n. The 
equilibrium conditions between bulk crystalline solid and liquid phases are imposed on the thermo¬ 
dynamic relationships for the PFC model to obtain a procedure for determining solid-liquid phase 
coexistence. The resulting procedure is found to be in agreement with the method commonly used 
in the PFC community, justifying the use of the common-tangent construction to determine solid- 
liquid phase coexistence in the PFC model. Finally, we apply the procedure to an eighth-order-fit 
(EOF) PFC model that has been parameterized to body-centered-cubic (bcc) Fe [Jaatinen et ah, 
Physical Review E 80, 031602 (2009)] to demonstrate the procedure as well as to develop physical 
intuition about the PFC input parameters. We demonstrate that the EOF-PFC model parameter¬ 
ization does not predict stable bcc structures with positive vacancy densities. This result suggests 
an alternative parameterization of the PFC model, which requires the primary peak position of 
the two-body direct correlation function to shift as a function of n. 


* vicchan@umich.edu 
i kpnirand@kmitl.ac.th 
i kthorn@umich.edu 


1 



I. INTRODUCTION 


The phase-field crystal (PFC) model is a simulation approach for studying phenomena 
that occur on atomic length and diffusive time scales. This is achieved by considering a 
free energy that is minimized by either a periodic order parameter profile, which represents 
a solid crystalline phase, or a constant order parameter profile, which represents a liquid 
phase nna. Such a formulation allows the PFC model to describe elastic and plastic 
deformation, multiple crystal orientations, and free surfaces in non-equilibrium processes 
[2J. Consequently, the model has been applied to investigate many important materials 
phenomena such as dislocation dynamics M, nucleation m, and grain-boundary-energy 
anisotropy Eli- 

The links between PFC model parameters to measurable quantities in experiments and 
atomistic simulations were made by Elder et al. |S] who showed that the PFC model can 
be derived from the classical density functional theory (cDFT) of freezing [TO], [IT]. This 
derivation provided a statistical mechanical interpretation of the PFC order parameter as 
an atomic-probability density, which is obtained by taking an ensemble average of the mi¬ 
croscopic particle density ra- The derivation also associated the bulk modulus and lattice 
spacing of a crystal to the curvature and position, respectively, of the first peak of the 
two-body direct correlation function (DCF), which can be obtained from experiments or 
atomistic simulations. 

Although the PFC model parameters have been linked to measurable quantities, the 
procedures for calculating equilibrium material properties from the PFC model are not 
straightforward ra because the thermodynamic interpretation of the PFC free energy has 
not been fully developed. In this paper, we present a thermodynamic interpretation for bulk 
phases of the PFC model. As a starting point, we follow the thermodynamic formalism 
for a crystalline system that was introduced by Larche and Cahn [QJ and was detailed 
in Voorhees and Johnson [13j to derive a thermodynamic relationship between the PFC 
free energy and thermodynamic state variables. We then apply the equilibrium conditions 
between a bulk crystal and liquid from Voorhees and Johnson ra to the thermodynamic 
relationship for the PFC model to obtain a thermodynamically consistent procedure for 
determining solid-liquid phase coexistence, which is demonstrated to be in agreement with 
the common-tangent construction commonly used in the PFC community [2l [T6l ITT| . Finally, 
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we apply this procedure to a PFC model parameterized for body-centered-cubic (bcc) Fe via 
an eighth-order fit (EOF) of the two-body DCF in Fourier space [IS]. The EOF-PFC model 
is used to demonstrate the procedure as well as to examine how the average value of the 
order parameter, n, and lattice spacing, a, are related to the number of atoms and vacancies 
in a crystal simulated by the PFC model. 

The paper is outlined as follows. In Section [TTJ we review the free energy of the PFC 
model. In Section [TTT] we describe the thermodynamics for bulk liquid and crystalline phases. 


In Section IV we derive the free-energy densities (FED) for the bulk liquid and crystalline 
phases, which serve as the basis for the thermodynamic interpretation of the PFC free 


energy. In Section VI we present the equilibrium conditions between a bulk crystal and 
liquid phase and apply these conditions to the thermodynamic relationship for the PFC 


model to obtain a procedure for determining solid-liquid phase coexistence. In Section VII 


this procedure is applied to the EOF-PFC model to demonstrate the procedure, as well as 
to develop physical intuition about the PFC model parameters. For this model, we derive 
an upper-bound expression for n, above which the vacancy density becomes negative. We 
further show that the EOF-PFC model does not stabilize bcc structures if values of n are 
below the upper bound. These results indicate a need for an alternative parameterization of 


the PFC model. Finally, in Section VIII we summarize the results of our work and present 
potential directions for future work. 


II. THE PFC FREE-ENERGY FUNCTIONAL 

The PFC model is based on a free-energy difference with respect to a reference liquid 
phase. The free energy is written in terms of an ideal-gas contribution, AJ-)d[n(r)], and an 
excess contribution, AJ ei [n(r)], 

AT[n(r)] = AJ r id [?r(r)] + A7' ca .[n(r)], (1) 

where A7),j[n(r)] is derived from a system of non-interacting particles and AJ r ea; [n(r)] con¬ 
tains the contribution from the interactions between particles. |9j The A symbols denote 
free-energy differences with respect to a reference liquid phase with (constant) atomic- 
probability density of p Q . The scaled dimensionless number density, n(r), is related to 
the atomic-probability density, p(r), by n(r) = p(r)/p 0 — 1. 
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Each term in Eq. (JTj) is written as an integral of a FED, 
AJ[n(r)] = [ A/ PFC (n(r))hr 

A/id(n(r)) + Af ex (n(r)) 


dr, 


( 2 ) 


where Af PFC (n(r)) is the FED of the PFC model, and Af id (n(r)) and Af ex (n(r )) are the 
ideal and excess contributions to Af PF c(n(r)), respectively. The subscripts V denote that 
the integrals are over the system volume. 

The ideal-gas FED, 


Af id (n(r)) = p 0 k B T 


n( r) 2 n(r) 3 n(r) 4 " 1 

—v- b t — -h 


6 


12 


(3) 


is obtained by approximating the Helmholtz free energy of an ideal gas with a polynomial 
expansion. The parameters at and bt are fitting constants, k B is the Boltzmann constant, 
and T is the temperature of the reference liquid phase. On the other hand, the excess FED, 

Af ex (n(r )) = - p ° k ^ T J n ( r )C (2) (| r - r |)ra(r')dr', (4) 


is obtained by expanding the excess Helmholtz free energy to include correlation functions 
up to the second order, i.e., the two-body DCF, C^ [12] • In writing the two-body DCF 
as C (2) (| r — r' |), an assumption has been made that the DCF is spherically symmetric.|9] 
Combining Eqs. (Jl| through ([I]), the PFC free energy is 

AJ zr [n(r)] = p 0 k B T J (a t - y*- bt ~^~ + ~^2 -^ j n ( r )C (2) G r - r' |)u(r')dr'^ dr. 

(5) 

AJ r [?r(r)] is minimized by n(r) that is equal to a constant value or that contains peaks with 
the periodicity of a crystalline lattice. Regions where n(r) is constant are considered to be 
in a liquid state, while those where n(r) is periodic are considered to be in a crystalline 
state. 

The calculation of Af ex (n(r)) is often performed in Fourier space, 

A/e*(n( r)) = - Po ^ BT n(r)^~i[C ,(2 )(fc)u(k)], (6) 

where the convolution theorem is used to efficiently evaluate the integral of Eq. Q as the 
inverse Fourier transform of the product of Fourier transforms. The notation ] denotes 
the inverse Fourier transform, k is a wave vector, k = |k|, and the hat symbols denote the 
Fourier transform of the respective quantities. 
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Thermodynamics describe the properties of systems that are in equilibrium. In the PFC 
model, an equilibrium density profile, n eq ( r), is obtained by relaxing n(r) via conserved 
dissipative dynamics PM, 


dn( r) = V 2^^K01 

dt 8n( r) 


(7) 


until a stationary state is reached. This state corresponds to the lowest energy state for the 
given constraints on the order-parameter average, n = (1/V) J n(r)dr, and lattice spacing, 
a, and will be referred to as the single-phase equilibrium state. In a single-phase equilibrium 
bulk liquid phase (i.e., away from any interfaces or boundaries), n eg (r) = ng“ lk,, (r) = n. 
Thus, a coarse-grained FED for the bulk liquid phase is A/ppjv* (ft) = Af id (n) + A f ex {n). 
On the other hand, in a single-phase equilibrium bulk crystalline phase, n eq { r) = n^ llk,c (r), 
which has a density profile that is periodic with a uniform amplitude for each value of n. 
Therefore, the free energy of a bulk crystalline phase corresponding to n^ q lk,c (r) is a function 
of n, a, and system volume, V. Consequently, the coarse-grained FED of a bulk crystalline 
phase is given by 


a />bulk,c 
PFC 


(n, a) 


lim 

V—>oo 


Ajr[ n bu lk , C ( r )] 

V 


( 8 ) 


where the limit indicates that the system volume is large enough such that the bulk phase 
is far away from any interfaces or boundaries. The definition in Eq. ([8]) shows that the PFC 
FED of a bulk crystalline phase is a function of n and a, and is defined in terms of V. Since 
n and a are coarse-grained variables, a system having n(r) = n^ lk,c (r) for the solid phase 
and n(r) = n k q lk ’ l (r) = n for the liquid phase can be described by the thermodynamics of 
bulk phases. 

The value of a that minimizes A/pp lk (n, o) for each n is denoted as a*. Its value is set 
by the position of the maximum of the primary (first) peak of the two-body DCF in Fourier 
space, k m , where a* oc k- 1 01]. Since the position of the primary peak in the two-body 
DCF is independent of n in the PFC model, the value of a* remains unchanged for all values 
of n. 

The PFC model has several variations that depend on the choice of the two-body DCF. 
A thorough review of different formulations and extensions of the PFC model is given in 
Ref. EE In this paper, we consider the PFC model proposed by Jaatinen et al. [18] because 
their model has been parameterized to an experimentally based two-body DCF that was 
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calculated from atomistic simulations of bcc Fe [21] . Their PFC model employs an ideal-gas 
FED with a t = 0.6917 and b t = 0.08540 and an EOF approximation of the two-body DCF, 
which has the form 

pocSfW = - r ’ (9) 

where 

r = E B =C m -C 0 -T , (10) 

o 

and k m , C 0 , C m , and C c are fitting constants that are defined in terms of the features 
of the two-body DCF in Fourier space: C 0 = PoC F q F ( 0), C m = p 0 C F Q F (k m ), and C c = 
d 2 (yPoC F Q F {k/dk 2 \k=k m (curvature at k = k m ). In Ref. Q21 the EOF-PFC model was pa¬ 
rameterized to the bulk modulus and solid-liquid interfacial energies of bcc Fe, where the fit¬ 
ting constants have the following values: k m = 2.985A -1 , C 0 = — 49A~ 3 , 1 — C m = 0.332A” 3 , 
C c = -10.40A- 1 , and p 0 = 0.080lA- 3 . 

III. THERMODYNAMICS FOR BULK PHASES 

In this section, we review the thermodynamics for bulk liquid and crystalline phases. 
First, we consider the thermodynamics for a bulk liquid phase, i.e., n(r) = ng” lk3 (r) = h, 
and provide definitions for the entropy, S l , hydrostatic pressure, P l , and chemical potential, 
p l A , for the bulk liquid phase, where the superscript l denotes quantities associated with 
the liquid phase. Next, we consider the thermodynamics of a bulk crystalline phase, i.e., 
?r(r) = 7ig“ lk,c (r), which incorporates the description of a lattice with sites that contain either 
atoms or vacancies, as introduced by Larche and Cahn [T^. The definitions for the entropy, 
S c , hydrostatic pressure, P c , and two different chemical potentials, p c A and p F , are provided 
for the bulk crystalline phase, where the superscript c denotes quantities associated with 
the crystalline phase. 

A. Bulk Liquid Phase 

The free energy of the PFC model is the Helmholtz free energy [9]. For the bulk liquid 
phase, the Helmholtz free energy is 

F l = F l {0\ V 1 , N l A ), (11) 
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where 0 l , V , and N l A are the temperature, volume, and the number of atoms in the bulk 
liquid, respectively. The differential form of F l is 

dF l = -S l d9 l - P l dV l + n l A dN l A (12) 


where 


S l = 


dF l 


de l 


p l = 


dF l 


V l ,N‘ 


dV l 


e l ,N‘ 


l l A — 


dF l 


dN\ 


(13) 


e l ,v l 


The integrated form of Eq. (12) is 


F\6 l , V\ N l A ) = -P l V + fi l A N l A 


(14) 


and the corresponding Gibbs-Duhem relation is 

0 = S l dd l - VdP 1 + N l A d/i l A . 


(15) 


Equations (14) and (15) will be used later to derive FEDs for the bulk liquid phase. 


B. Bulk Crystalline Phase 


For a bulk crystalline solid, a lattice is employed to represent the spatially periodic 
structure of a crystal, where each lattice site contains either an atom or a vacancy 01 
For a one-component bulk crystalline phase, the total number of lattice sites, 1V£, is related 
to the number of atoms, N A , and vacancies, Ny, by 

N c l = N c a + N£. (16) 


As discussed in Voorhees and Johnson ra. any two of the three thermodynamic variables 


in Eq. (16) can be used to describe the thermodynamic state of a one-component crystal. 
In this work we consider the Helmholtz free energy as a function of N A and N £, 


F c 


F%e\V\N%,NQ, 


(17) 


where 9 C and V c are the temperature and the crystal volume, respectively. The differential 
form of F c is 

dF c = —S c d9 c - P c dV c + fj, c A dN c A + dldNl (18) 
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where 


S c = 


OF' 


dO c 


d F c 


dF c 


c 9 F c 

pc — 

’ dV c 

V c ,N\,Nl uv 

1 

e c ,N^,Ni 

III 

1 

9 c ,V c ,Nl 

^ L ~ dN c L 


e c ,v c ,m 


(19) 


The partial derivative that defines n A provides the energy change due to the addition or 
removal of an atom while the number of lattice sites and crystal volume are held constant 
for an isothermal system. The definition of fi c A requires a vacancy to be eliminated when an 
atom is added to the crystal and a vacancy to be generated when an atom is removed from 
the crystal. This chemical potential is hereafter referred to as the diffusion potential ra 
to distinguish from the chemical potential of atoms for the bulk liquid phase. On the other 
hand, the partial derivative that defines /j, c l provides the energy change due to a change 
in the number of lattice sites while the crystal volume and the number of atoms are held 
constant for an isothermal system. A lattice site can be added by moving an atom within 
the crystal to the surface while simultaneously creating a vacancy within the crystal. This 
process will cause an increase in the pressure when the crystal volume is held constant and 
the partial molar volume of the vacancy is nonzero [15] . 


The integrated form of Eq. (18) is 


F C (9 C , V c , N%, Ni) = —P C V C + fi c A N c A + ^ C L N C L 
and the corresponding Gibbs-Duhem relation is 

0 = S c d9 c - V c dP c + N c A d^ c A + Nidul, 


( 20 ) 


( 21 ) 


the derivation of which is presented in appendix |A| Equations (20) and (21) will be used 
later to derive FEDs for a bulk crystalline phase. 


IV. FREE-ENERGY DENSITIES FOR BULK PHASES 

In this section, we derive FEDs that are defined on a reference (undeformed) volume, 
V', and the (potentially deformed) system volumes, V 1 and V c , for the bulk liquid and 
crystalline phases, respectively. FEDs defined on V are referred to as reference-volume 
FEDs, while those defined on V 1 or V c are referred to as system-volume FEDs. By defining 
the reference-volume FED for the bulk crystalline phase, we make an important distinction 
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between two sources of pressure change: mechanical and configurational forces [15]. The 
system-volume FEDs are used to develop a thermodynamic relationship between the PFC 
FED and thermodynamic state variables. 


A. Bulk Liquid Phase 


The reference-volume FED for the bulk liquid phase, denoted as /(,,, is obtained by 
dividing F l by V', 

T?l 

( 22 ) 


fi'=y7= ~P‘J‘ + 


where 


,1 V ‘ fl n ‘a 
J = V’ 


(23) 


The variable J l describes volume change due to hydrostatic pressure and p A is the atomic 
density of the liquid phase defined on V', as denoted by the prime symbol. Similarly, the 


Gibbs-Duhem relation in Eq. (15) is divided by V to obtain 


0 = s l v ,d6 l - J l dP l + p n A dp l A , 


(24) 


where s l v , = S l /V. Differentiating Eq. (22) and subtracting Eq. (24) gives an expression for 

dfv>, 

dfy, = —s l v ,d6 l - PdJ 1 + p l A dpf{, (25) 

where /y, is a function of natural variables 9\ J l , and p A . 

Alternatively, the system-volume FED for the bulk liquid phase, as denoted by /y, is 
obtained by dividing F l by lA, 


fv = — -P l + d l AP l Af 


(26) 


where p l A = N l A /V l . Similarly, the Gibbs-Duhem relation in Eq. (15) is divided by V 1 to 
obtain 

0 = s l v dO l - dP l + p\dp\, (27) 

where s l v = S l /V l . The variables s l v and p l A are dehned on the system volume of the bulk 


liquid phase, V . Differentiating Eq. (26) and subtracting Eq. (27) gives a relationship for 

dfl 

dfv = —s l v d6 l + p l A dp l A , (28) 
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where /y is a function of natural variables 9 l and p l A . Equation (28) is defined in terms of 
the system volume of the liquid phase and is a function of p l A , which is related to n of the 


PFC FED for a bulk liquid via n = p l A /p 0 — 1. Therefore, Eq. (28) will be employed to 
develop a thermodynamic interpretation of the PFC free energy for the bulk liquid phase. 


B. Bulk Crystalline Phase 


The reference-volume FED for the bulk crystalline phase, denoted as /■£,, is determined 
by dividing F c with V', 


where 


fv'=y- -P C J C + PaPa + pIpI 


JC _ V c /c _ N% , c _ N c l 

- V' ’ Pa ~ V ’ Pl ~~ V' ' 


(29) 


(30) 


Similar to the bulk liquid phase, J c describes volume change due to hydrostatic pressure in 
a bulk crystalline phase, and p A and p'f are the atomic and lattice densities, respectively, 


that are defined on the reference volume. The Gibbs-Duhem relation in Eq. (21) is divided 
by V to obtain 

0 = s c v ,d6 c - J c dP c + p ( A dp c A + p'ldn c L , (31) 


where Sy, = S c /V. Differentiating Eq. (29) and subtracting Eq. (31) gives an expression for 
dfv o 

dfc, = - s c v d0 c - P c dJ c + p c A dp'Z + p c L dp'l, (32) 


where /y, is a function of natural variables 9 C , J c , p A , and p'[. 

Changing the volume of an isothermal bulk crystalline system (reflected by a change in 
J c because V is constant) while keeping the mass (equivalent to fixing p A ) and number 
of lattice sites (equivalent to fixing p'l) constant will cause the pressure to change. This 
pressure change arises from deforming the system with a mechanical force. Alternatively, 
as mentioned earlier, a pressure change also arises when the number of lattice sites change 
(reflected by a change in p'f) while the crystal volume and mass are constant. To understand 
this latter type of pressure change, consider a thought experiment where a crystal is enclosed 
by a rigid wall that does not allow mass transfer (fixed V c and p' A ). When a lattice site 
(vacancy) is added to the system, the constraint imposed by the walls prevent a volume 
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change, and thus results in a pressure change. This type of pressure change arises from a 
configurational force [T5]. Therefore, the thermodynamic framework for crystalline solids 
described above allows us to distinguish between pressure changes due to mechanical and 
configurational forces m 

As mentioned earlier, the PFC FED is defined on the system volume. Therefore, a 
system-volume FED for a bulk crystalline phase will be most appropriate for developing 
a thermodynamic relationship for the PFC model of a bulk crystal. To reformulate the 


reference-volume FED, Eq. (29), and Eq. (31) in terms of the system volume, we divide 


both sides of the equations by J c (see Eq. (30) for definition): 


ff c - P c , c Pa , c Pl 

ff v = JP = ~ P +Paj-c + Plj~c 


and 


0 = -^d0 c 
J c 


dP c + f c dA + P fdp\, 


(33) 


(34) 


where ffy is one expression for the system-volume FED. The division by J c above maps 
p A , Pl, and Sy, to their system-volume counterparts: p c A , p c L , and Sy, respectively. The 


differential form of Eq. (33) is obtained by taking its derivative and subtracting Eq. (34), 


dfi = - a fc d0C + Jd dp '° A + jPti + ^ + Ap '° l) 


(35) 


where it can now be observed that ffy is a function of natural variables 9 C , p A , p'f, and 1/J C . 
Therefore, the FED denoted by ffy is a system-volume FED that is a function of densities 
defined on the reference volume. 

In order to obtain a FED with independent variables that match those of the PFC FED, 
the variables p A and p'f are related to p c A and p c L by the chain rule 

df c A = d(p A J c ) = J c dp c A + p c A dJ c and dp'l = d(p c L J c ) = J c dp c L + p c L dJ c . (36) 


Substituting Eq. (36) into Eq. (35), one obtains 


dfi = -S^d6 c + ir A dp c A + fi c L dp c L) 


(37) 


where we have used the relationships 




dJ c = -(J c ) 2 d(l/J c ) and 4 = 


(38) 
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The same expression can be obtained by starting from a FED based on the system volume; 
however, we derive it in this manner to explicitly illustrate the connection between reference- 
volume and system-volume variables. This expression of the system-volume FED for the 


bulk crystal, /y, given in Eq. (37) is now a function of natural variables 9 C , p A , and p c L . By 
expressing p c L in terms of the volume of a lattice site, V£, 

1 


Pl = 


V£’ 


(39) 


Eq. (37) is related to the lattice spacing via the unit-cell volume, V£, which is written in 
terms of V£ using the number of lattice sites per unit cell, x° = T£/V£. The value of \ c 
depends on the lattice structure. For example, a bcc structure, which contains 2 lattice sites 
per unit cell, has y c = 2. 


Equation (39) is used to express d/£ as 


rf/'v = -s c v d6 c + II \df) c A + x c dl d ( l/ v c). 


(40) 


where /£ is a function of natural variables 9 C , p A , and 1/V£. Since is a sole function of 
a (e.g., for a bcc structure, = a 3 ), the FED representation in Eq. (40) is a function of 


natural variables that correspond to those of the PFC FED of the bulk crystal. Therefore, 


Eq. (40) is employed in developing a thermodynamic interpretation of the PFC free energy 


for a bulk crystalline phase. The integrated form of Eq. (|40|) is derived from Eqs. (33) and 

(41) 


(39) to be 


fv = - pC + PaPa + 


X°P C L 


V, 


c 


and Eqs. (36) and (39) are substituted into Eq. (34) to obtain 


0 = - s c v d9 c - dP c + p c A dpf A + ^-dp c L . 

v c 


(42) 


The expressions in Eqs. (40) and (41) form the basis for deriving a thermodynamic relation¬ 


ship for a bulk crystalline phase in the PFC model. 


As observed from Eq. (40), p A = d /£/ dp A \e c ,v£ an d thus changing p A while holding 
and 9 C constant is equivalent to changing the number of atoms while the number of lattice 


sites and crystal volume are held constant (see Eq. (19)). A negative value of p A indicates the 
presence of a driving force for adding an atom into a vacant lattice site. On the other hand, 
a positive value of p, A indicates the presence of a driving force for removing an atom from an 
occupied lattice site. Additionally, it can be observed that x°Pl = ^fv/^(^/^c)W,p c A an d 
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thus changing V<h (equivalent to changing lattice spacing) while holding p A and 9 C constant 
is equivalent to changing the number of lattice sites while the crystal volume and mass are 
held constant. Therefore, the process of changing Vq while holding p c A constant gives rise 
to a configurational force, just as in the case for changing the number of lattice sites as 


indicated in Eq. (19). This point will be further examined in the next section after the FED 


in Eq. (40) is linked to the PFC FED. 


V. ISOTHERMAL THERMODYNAMIC RELATIONSHIPS FOR THE PHASE- 
FIELD CRYSTAL MODEL 

In this section, we use the system-volume FEDs derived in the previous section to develop 
the relationship of the PFC FED to thermodynamic state variables for the bulk liquid and 
crystalline phases. Here we limit the scope to isothermal systems, which correspond to a 
fixed two-body DCF in the PFC model. Therefore, d6 = 0 for our formulation below. 


A. Bulk Liquid Phase 


Since the PFC model is based on a free-energy difference from a reference liquid phase, 


the PFC FED of a bulk liquid phase is related to /{, in Eq. (|26j) by 

A/£*'(«) = /{, - /„, 


(43) 


where / 0 is the Helmholtz FED of the reference liquid phase with a density of p 0 . Note that 
/o remains constant as n changes because the two-body DCF is taken at the reference state 


and is assumed to be independent of n. Furthermore, n is related to p A of Eq. (26) by 


p A = np 0 + po- 


Equations (43) and (44) are combined with Eq. (26) to obtain 


^■fpFc\ n ) = fv ~ fo — —P l + PoPa ( n + 1) — /o- 


(44) 


(45) 


Similarly, Eq. (44) is combined with Eq. (27) to obtain 

0 = -dP l + p 0 (h + l)dp l A . (46) 

Differentiating Eq. (45) and subtracting Eq. pbj ) gives an expression for d(A/ppjV(h)), 

d( A /pFC ''(w)) = d (/y — f 0 ) = p 0 p l A dn. (47) 
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B. Bulk Crystalline Phase 


The PFC FED of a bulk crystalline phase is related to /y in Eq. (41) by 


A/S''(n,a) = /5-/„, 


( 48 ) 


and n is related to p A of Eq. (41) by 


Pa = n Po + Po- 


(49) 


Equations (48) and (49) are combined with Eq. (41) to obtain 


A/ P “ k,C («. a) = fi-f 0 = -P‘ + p„ii c A (n + 1) + 


X Pl 
V; 


/o- 


(50) 


c 


Similarly, Eq. (49) is combined with Eq. (42) to obtain 


X 


0 = — dP c + p 0 (n + 1 )dp A + —dp c L . 


V, 


(51) 


c 


Differentiating Eq. (50) and subtracting Eq. (51) gives an expression for d(Afpp^’ c (n, a)), 
d(A/“ 0 k ' c (fi, a)) = iUi- fa) = PaVf A dn + *%<i( 1/V&), (52) 


where is a function of a. 


Equation (52) is central to this work because it links the PFC FED to p A and p,y, which 
are chemical potentials that correspond to different thermodynamic processes. As seen in 


Eq. (52), p c A = 1/p 0 (<9A/ppQ ,c / dn\ V c) and therefore varying n when the temperature and 


lattice spacing are held constant in the PFC model is equivalent to changing the number 
of atoms while 9 C , V c , and 1V£ are fixed. Similarly, since p c L = l/x c (0A/pp 1 ^’ c /c}(l/V^)|f i ), 
varying lattice spacing when temperature and n are held constant in the PFC model is 
equivalent to changing the number of lattice sites while 9 C , V c , and N A are fixed. 

As discussed earlier, the addition or removal of lattice sites while the crystal volume is 
held constant gives rise to a pressure change due to a configurational force. This point 
has not previously been elucidated, and the resulting pressure has been instead attributed 
to pressure change due to mechanical forces, leading to an improper procedure for elastic 
constant calculations with the PFC model [2 DU [23]. A thermodynamically consistent 
procedure for calculating elastic constants was developed in our previous work [13] , and the 
framework above further validates our approach. 
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VI. SOLID-LIQUID PHASE COEXISTENCE IN THE PFC MODEL 


In this section, we apply the thermodynamic relationship developed in the previous sec¬ 
tion to derive a procedure for determining solid-liquid phase coexistence in the PFC model. 
First, we describe the equilibrium conditions presented in Voorhees and Johnson [15j be¬ 
tween bulk crystalline and liquid phases, which involve constraints on P, p,A, and p c L . These 
equilibrium conditions are then imposed on the thermodynamic relationship for the PFC 


model, Eqs. (47) and (52), to obtain a thermodynamically consistent procedure for deter¬ 
mining solid-liquid phase coexistence. 


A. Solid-Liquid Phase Coexistence in the Phase-Field Crystal Model 


As described in Voorhees and Johnson [TS], a bulk crystalline and liquid phase are in 
equilibrium when the variation in the total energy of a system vanishes. They showed 
that the variation in the system energy vanishes when thermal, mechanical, and chemical 
equilibria are achieved. For a one-component bulk crystalline phase in equilibrium with a 
bulk liquid phase, the equilibrium conditions are 

0 C = 0\ P c = P l , /4 = /4 »L = 0- (53) 


The last condition of p c L = 0 is unique to crystalline solids and indicates that there is no 
driving force for an addition or removal of a lattice site. This condition can be used to 
calculate equilibrium vacancy concentration for a crystal in equilibrium with a liquid phase 

m- 


The condition for thermal equilibrium is met when considering an isothermal system. 


Furthermore, p c L = 0 when the FED of the bulk crystalline phase, A/ppQ ,c (n, a), is relaxed 


with respect to the lattice spacing while n and 6 C are held constant (see Eq. (52)). Therefore, 
the task of finding solid-liquid coexistence lies in satisfying the equilibrium conditions, P c = 
P l and p c A = p l A . 

A relationship for pressure is obtained by rearranging the expressions for the PFC FEDs 


in Eqs. (45) and (50): 


-P l 


A/ppc ' (ra) - polJ- l A {n T 1) + /o and - P c 


A /pfc ' X a)-pop A (n + 1) 




F/o, 

(54) 
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respectively, where the determination of P l and P c require a knowledge of /o. However, 
since both the liquid and solid phases have the same reference state, / 0 will cancel when 
equating their pressures. By using 


r\ \ pbulk.Z 

, _ dA f PFC : 
Pop a ~ 


dn 


(55) 


obtained from Eq. (47), P l in terms of A/p F( V(n) and n is expressed as 


bulk,/ 


- P‘ = A/p^'tn) - Ahic 

= y l -s l + f 0 


(n + 1) + /o 


e l 


(56) 


where 


l _ /\ pbulk,;, 

PFC 1 


y = a/] 


n) 


r\ A r*bulk,/ 

cl - dA/ppc’ 


S l n and S = 


dn 


e l 


are the y-intercept and slope of a A/ppj’y^n) vs. n curve, respectively, when 6 l is constant. 
A schematic of Eq. (57) is shown in Fig. [I] 


(57) 



FIG. 1. Schematic of Eq. (57) on a liquid FED curve (dashed line) with tangent line 


at the point n = n *, which is marked with “x”. The schematic for Eq. (60) is similar, 
FED curve of a bulk crystalline phase. 


(solid line) 
but for the 


Similarly, by using 


r\ A /?bulk,C 

_ dA / PFC ’ 
PoPa ~ 


dn 




( 58 ) 
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obtained from Eq. (52), P c in terms of A/pp^’ c (n, a) and n is expressed as 


P ' — A/pfc 

= y°-S 


bulk.c/— \ 

l n , a) 


dAf, 


bulk,c 

PFC 




V, 


dn 


+ fo, 


[n + 1) 


X c /4 


0 C ,» C L 


V, 


+ /o 


c 


c 


where 


^ = A/j 


bulk,c/ 

PFC 


n, a 


S c n, and S c = 


dAf, 


bulk,c 

PFC 


dn 


( 59 ) 


(60) 


»',et 


are the y-intercept and slope of a A/pp^’ c vs. n curve, respectively, when fi c L and 9 C are 


constant. Note that the condition /i c L = 0 has not been imposed in Eq. (59), but will be 


applied later to fulfill the equilibrium conditions listed in Eq. (53). 


According to Eqs. (56) and (59), mechanical equilibrium for solid-liquid phase coexis¬ 
tence is fulfilled when y c = y l , S c = S l , and )i c L = 0. Furthermore, the condition, /i c A = /J A , 
for chemical equilibrium is also satisfied when S c = S l . Therefore, satisfying the common- 
tangent condition (i.e., the y-intercept and slope of the two curves are equal) for the A/pp 1 ^’ 0 
vs. n and fppQ l (n) vs. n curves with fi c L = 0 fulfills the equilibrium conditions for solid-liquid 
phase coexistence. This procedure, which involves thermal, mechanical, and chemical equi¬ 
libria, is in agreement with the common-tangent construction commonly used for calculating 
solid-liquid phase coexistence in the PFC model [2j 16], [IT]. Therefore, the above analysis 
justifies the use of the common-tangent construction to determine solid-liquid phase coexis¬ 
tence in the PFC model. 

As noted earlier, the condition of /ip = 0 is accomplished in the PFC model by minimizing 
the PFC FED for the bulk crystalline phase with respect to Vq while 6 C and ft are held 


constant, as seen in Eq. (52). Therefore, the minimization of the density profile according to 
the PFC FED with respect to lattice spacing for a crystalline solid, which is conventionally 
done in the PFC model [16] , is a necessary step for obtaining the state in which fi c L = 0 and 
for determining solid-liquid phase coexistence. 


VII. APPLICATION TO EOF-PFC MODEL 

In this section, we apply the thermodynamic relationships for the PFC model developed 
in the previous section to the EOF-PFC model (see Section |nj) to demonstrate the procedure 
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for determining phase coexistence, as well as to develop physical intuition about the PFC 
input parameters. 


A. Free-Energy Density Curves and Phase Coexistence for bcc Fe 

We calculate the solid and liquid FED curves for the EOF-PFC model with the fitting 
parameters of Ref. [18j which are listed in Section |TT} The dimensionless PFC FED for the 
solid and liquid phases are plotted as functions of n in Fig. [2] The equilibrium FED for 


e 

II 

e 

s 

<1 



n 


FIG. 2. Plot of the dimensionless PFC FEDs of the EOF-PFC model for the solid (solid red 
line) and liquid (dashed blue line) phases as a function of n. Each point on the solid FED curve is 
minimized with respect to lattice spacing, and thus satisfies /.i c L = 0. The dimensionless coexistence 
number density for the solid, n s , and liquid, n/, phases are marked with “x” marks. 

each value of n was calculated by relaxing a one-mode approximation of bcc Fe according 
to the PFC FED via Eq. ([?]), as done in Ref. [TS] Each point in the solid FED curve of Fig. 
[2] satisfies fi c L = 0, which is achieved by minimizing the PFC FED with respect to lattice 
spacing. As described in Section [Til the value of a that minimizes the PFC FED, denoted as 
a*, remains constant because the position of the primary peak in the two-body DCF of the 
PFC model is assumed to be independent of n. For the calculations in Fig. [ 2 J a* = 2.978A 
for all values of ft. 

The scaled dimensionless coexistence number densities for the solid, n s , and liquid, h;, 
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phases were determined with a common-tangent construction on the FED curves in Fig. 
[ 2 } The values for solid and liquid coexistence densities are n s = 9.17 x 1CT 3 and ni = 
—2.49 x ICE 2 , respectively. These values are in agreement with those presented in Ref. [IB) 


B. Diffusion Potential 


Figure |3] shows the diffusion potential calculated from Eq. (58) for the single-phase solid 
region of the EOF-PFC model, n > n s . As can be observed, fi A transitions from a negative 



FIG. 3. Plot of diffusion potential, n A , for n > n s , where only solid is stable (see Fig. [ 2 ]). The 
dashed vertical line corresponds to n = 3.54 x 10“ 2 , which is the value of n where /i A = 0 (denoted 
by horizontal dashed line). 

to a positive value at n = 3.54 x 10” 2 . As described earlier, n c A represents the energy change 
due to the addition or removal of atoms in an isothermal system when the number of lattice 
sites and the crystal volume are held constant. When /j, a < 0 (left of dashed vertical line in 
Fig. |3[. there is a driving force for adding an atom into a vacant lattice site, which decreases 
as n increases. On the other hand, when fj, A > 0 (right of dashed vertical line in Fig. [3]), 
there is a driving force for removing an atom from an occupied lattice site, which decreases 
as n decreases. When fi A — 0, there is no driving force for adding or removing atoms to 
and from lattice sites, and the condition fi c L = 0 allows us to determine the single-phase 
equilibrium vacancy density. 
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In the solid-liquid coexistence region (fq < n < n s ), chemical equilibrium requires that 
the chemical potentials of the two phases are equal, which will give rise to nonzero diffusion 
potentials. Thus, the equilibrium vacancy density is determined from p c A = jJ A and p c L = 0. 


C. An Upper Bound for n 


As discussed earlier, the number of atoms, vacancies, and lattice sites of a crystal are 


related to each other by Eq. (16). As a result, the vacancy density, py, can be expressed in 


terms of p A and p c L as py = p c L — p A . By substituting Eqs. (39) and (49), py is expressed in 


terms of u and Vq as 


Pv = ~ Po{n + l). 

v c 


(61) 


Since Vq is constant for all n values (see Section VII A), an upper bound for n arises when 
the value py is specified. An upper bound for u, h max , is obtained when Py = 0 (i.e., crystal 
with no vacancies), 

hmax = - I- (62) 

V c c pQ 

For n > U max , the vacancy density takes a negative value, which is unphysical. 

The upper bound for the EOF-PFC model with the fitting parameters described in Section 
HH is u m ax = —5.46 x 10 —2 , where y c = 2 for a bcc structure, Vq = (2.978A) 3 , and p 0 — 
0.0801A -3 . p8] Surprisingly, U max < ni, where the solid phase is unstable. Therefore, the 
EOF-PFC model parameterized in Ref. [18] does not stabilize a bcc Fe structure with py > 0. 

A potential interpretation of py < 0 is the presence of mobile interstitials. However, 
further investigation is needed to examine this possibility and its validity. In this work, 
we simply consider this case as an artifact of the model parameterization and proceed to 
suggest potential solutions. For example, changing the correlation function, as well as the 
parameterization of A/^(n(r)) in Eq. ([3]), can change the stability of the solid phase such 
that u; < U m ax. Another approach is to require the position of the primary peak of the 
two-body DCF, k m , to be a function of u, 

k m (n ) = 2ttV 1 2 + nri 2 + u 2 (V£.)“ 5 


where 


= 27ta/Z 2 + m 2 + u 2 


V c c (n) = 


X 


[p c v + Po(n+ 1 )] 


r 


\p c v + Po(n+ 1 )] 


(63) 

(64) 
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is obtained by rearranging Eq. (61), and l, m , and n are the Miller indices of the primary 
family of planes (e.g., I = 1, m = 1, and n — 0 for the bcc structure). In this case, the upper 


bound in n no longer arises. Equation (63) also allows the direct control of py, which has 
not been previously possible. It also indicates that the manner in which the two-body DCF 
changes with n depends on the crystal structure via x° and the Miller indices. Equation 


(63) should only be applied when the change in n is due to the addition or removal of an 
atom while 9 C , V c , and 7V£ are held constant. This corresponds to changing n while holding 
9 C and constant in the PFC model. Note that, in the case where the change in n is due 
to a change in V c while N\ is held constant, the value of k m must remain fixed in order to 
apply a mechanical force, which gives rise to a pressure change. 

We point out that a change in k m reflects a change in the liquid reference state. Therefore, 


a parameterization that requires k m to be a function of n (as in Eq. (63)) requires po and 
the liquid reference pressure, P°, to change with n. Since ft is a function of po, a relationship 
for k m as a function n can only be obtained from the dependence of k m on p 0 , which must 
be determined from atomistic simulations. 


VIII. SUMMARY AND DISCUSSION 

In this paper we have applied the thermodynamic formalism for crystalline solids of Larche 
and Cahn [Q] to develop a thermodynamic relationship between the PFC free energy and 
thermodynamic state variables. This relationship allows us to examine the thermodynamic 
processes associated with varying the PFC model parameters. We showed that varying 
n while keeping the unit-cell volume, (and thus lattice spacing, a), and temperature, 
9 C , of the bulk crystalline phase constant in a PFC simulation reflects the thermodynamic 
process of adding or removing atoms to and from lattice sites. Furthermore, changing 
the computational size of a PFC simulation, while keeping n and 9 C constant, reflects the 
thermodynamic process of adding or removing lattice sites. 

The equilibrium conditions between bulk crystalline solid and liquid phases were then 
imposed on the thermodynamic relationships for the PFC model to obtain a procedure for 
determining solid-liquid phase coexistence, which we found to be in agreement with the 
method commonly used in the PFC community. By using the procedure, we found that no 
stable bcc phase with a vacancy density greater than or equal to zero exists for the EOF-PFC 
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model that has been parameterized to bcc Fe [18]. Therefore, we proposed an alternative 
parameterization of the EOF-PFC model, which requires the position of the primary peak 
of the two-body DCF to be a function of n. The implementation of this parameterization 
will be a topic of future work. 

Although several thermodynamic processes associated with changing PFC input parame¬ 
ters (i.e., n and a) were elucidated for bulk phases in this work, a thermodynamic framework 
for interfaces in the PFC model remains to be developed. An extension of this work that 
considers interfaces is needed in order to gain a rigorous, quantitative understanding of 
PFC simulation results that contain interfaces and grain boundaries. Furthermore, since we 
do not consider the equilibrium conditions between different crystal phases, the use of the 
common-tangent construction to determine solid-solid phase coexistence ESI HZ! remains 
to be verified. The extension of this work to systems containing interfaces and different 
crystal phases will be the topics of future investigations. 
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Appendix A: Gibbs-Duhem Relation for Bulk Crystalline Solid 

We derive the Gibbs-Duhem relation for a bulk crystalline solid following Voorhees and 
Johnson H5j. The internal energy of a one-component crystal is a function of S c , V c , 
and iV£, 

E c = E c (S c ,V c ,N c a ,N c l ). (Al) 

The differential form of E c is 

dE c = 9 c dS - P c dV c + fj, c A dN c A + n c r dNl, (A2) 
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where 


e c = 


8E C 


8S C 


nc dE c 


c 9E C 


8E C 

P — 

’ ~ 8V C 

V c ,N%,Nl uy 

S c ,N%,Nl 

- 8N% 

5 

5 C ,V C ,./V£ 

III 

a. 


syv c , m 

(A3) 


Since Eq. (Al) is a homogenous function of degree one, Eq. (A2) yields 


E C (S C , V c , N c a , N c l ) = S C 9 C - P C V C + ijl c a N c a + fi c L N c L . 


(A4) 


Differentiating Eq. (A4) and subtracting from Eq. (A2) gives us the Gibbs-Duhem relation 


for a bulk crystalline solid, 


0 = S c dO c - V c dP c + N%dfi c A + N c L dn c L . 


(AS) 
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